Agent-based SLIR Modelling of Powdery Mildew Infection Spread in Grapevines

This notebook presents an agent-based model in combination with a SIR-like compartmental model of powdery mildew infection spread among grapevines in a spatial grid.

This is a study extending on a simple exploration done by Paul Petersik, with added details on consideration of infection and recovery transition and different aspects reflecting spatial heterogeneity displayed in the process of fungal infection spread.

Modelling and visualisation is done in python with AgentPy.

Preparation

AgentPy provides a basic framework for spatial modelling and a simple-to-implement visualisation solution. It is used in this study to establish a cellular automaton that simulates the changes of state of grapevines in the vineyard.

A popular alternative is Mesa.

# Model design
import agentpy as ap

# Visualization
import matplotlib.pyplot as plt 
import IPython

# General Utility
import random

Model Definition

This section details the models used in this study and a guide to walk through the python implementation of the model.

Model compartments

This study divides the grapevine population into 4 compartment:

1. Susceptible

   • Susceptible agents represent health grapevine individuals that are susceptible to potential fungal infection. Most agents in the simulation start in the susceptible compartment and transition into other compartments as powdery mildew spreads.

2. Latent

   • Latent agents are infected individuals without visible signs. Agents of this compartment are infected with powdery mildew but the infectants have not developed the capability of secondary spore production. During this process, the agent is dormant and progresses gradually to become actively infectious in the next state.

3. Infected

   • Agents fall into the infected compartment once they exit the latent period. The grapevine becomes infectious to its surrounding vines, while the fungi infection develops further in the host.

4. Recovered

   • Agents of the recovered compartment are grapevine individuals that have been previously infected and have recovered from the infection via either natural or manual treatment. These individuals gain higher resistance against infection from their recovery.

A class AgentState is created to host the 4 states/compartments for the model.

class AgentState():
    SUSCEPTIBLE = 0
    LATENT = 1
    INFECTED = 2
    RECOVERED = 3

Cellular automaton

The SLIR (Susceptible-Latent-Infected-Recovered) model runs inside a cellular automaton, with a collection of agents spatially represented as a grid. In context, the grid represents a sample field of grape vines and each cell, hosting an agent, represents an individual grape vine. As the automaton runs and iterates over generations, the agents updates their status according the pass of time and their surrounding circumstances.

The automaton has two main tasks: initialisation and stepping through generations.

Initialisation

The setup function is run when the model/automaton gets initialised. It prepares the grid and its agents, before the model starts to run, with an initial generation. It starts by creating the specified number of agents in this model and assigns them to the grid, filling on possible cells. The agents are then prepared with corresponding initial attributes, differing according to their spacial distribution. Such differences reflects spacial heterogeneity in grape vines' growth condition and inert ability.

The grapes vines are initialised as either part of the susceptible population or the resistant population. Such difference will be explained in the next section with the logic where spacial heterogeneity is introduced. However, both these compartments represent healthy grape vines. In order to set the automaton running, some infected individuals are randomly planted into the grid. These infected vines are in latent state so that their symptom will develop gradually until active and start to infect neighbouring vines, spreading the disease.

Stepping through generations

The model moves forward in time through generations of agents by calling the step function. The method is called in every iteration of generation, reflecting changes happening to grapevines in the vineyard. The function determines the next state of each individual agents in the grid, by first finding its neighbourhood. The state of a grapevine is influenced by infectious state of its surrounding vines in the spreading process of powdery mildew infection. The specific state-transition flow is detailed in the Agent section.

Once the next state of each agent is determined, the model sets the next state of each agent as its current state, signifying that a generation has passed and the necessary state transitions happening during this period have been done.

The model stops before reaching the maximum number of iterations when there are no more latent or infected agents in the vineyard. Since without these agents as active players, all agents in the grid remains static and simulation virtually ceases, no matter how time passes.

model_generation is a global variable that is updated by the automaton and tracks the current timestamp of the simulation. This variable is used later for time-related agent state transitions.

model_generation = 0

class VineInfectionModel(ap.Model):

    def setup(self):
        # Create agents
        vines = self.agents = ap.AgentList(self, self.p.size**2)

        # Create grid and add agents to the grid
        self.vineyard = ap.Grid(self, (self.p.size, self.p.size))
        self.vineyard.add_agents(vines)

        # Initiate a state variable for all agents
        set_up_agent_variable_attributes(self.vineyard)

        # Initialise model with infected vines
        initial_infected = self.vineyard.agents.to_list().random(
            self.p.initial_infection_count)
        initial_infected.state = AgentState.LATENT


    def step(self):
        global model_generation

        # Calculate the next state of each agent
        for vine in self.agents:
            neighboring_vines = self.vineyard.neighbors(vine).to_list()
            update_agent(vine, neighboring_vines)

        # Update the state of all agents
        for vine in self.agents:
            vine.state = vine.state_next

        # Update model generation count
        model_generation = self.t
        
        # Stop simulation if no active agent is left
        latent_agent_count = len(self.agents.select(self.agents.state == AgentState.LATENT))
        infected_agent_count = len(self.agents.select(self.agents.state == AgentState.INFECTED))
        if (latent_agent_count + infected_agent_count) == 0:
            self.stop()

Introducing spatial heterogeneity

Spatial heterogeneity is a main area of focus of this study. Measures are introduced to this simulation to simulate and reflect how the spread patterns of infection is affected spatially due to various factors.

In the current implementation, some of the main influencer of spread are:

• Intercropped resistant cultivars
• Variable infection latency period
• ...

Intercropping setup

The entire grid of vinegrapes are intercropped with a cultivar with more resistance to catch infection. This is done in order to setup barriers in between groups of normal grapevines to break the wave of fungal spread and slow down disease propagation.

The resistant cultivar perform much better against the infectant. However, it is not realistic to have too many such individuals in the overall population, as mentioned by Paul Petersik.

So, you may now ask the question: “Why not only plant resistant varieties?” The problem is, most consumers do not really like the wine made from resistant varieties (yet).

The following function demonstrates how intercropping is implemented in this model. The function uses the coordinates of the agent to determine whether the agent of interest is a normal grapevine cultivar or an intercropped resistant cultivar. Note that the pattern of intercropping is a evenly distributed row-by-row setup. As a result, whether an agent is intercropped resistant cultivar can be deduced by taking the remainder of the total pattern distance.

The above image visualises the intercropping setup, where light green spaces represent the normal grapevine cultivar and the darker green spaces represent the resistant cultivar.

intercropping_distance = 5
intercropping_row_count = 5

def is_intercropped_resistant_cultivar(position):
    x_coordinate = position[0]
    y_coordinate = position[1]
    return x_coordinate % (intercropping_distance + intercropping_row_count) >= intercropping_distance

Spacial variables

Latency duration: a gaussian distribution

Infectivity offset: cultivar based variable

default_infectivity_offset = 10
resistant_infectivity_offset = 180

default_latency_duration = 8

def set_up_agent_variable_attributes(grid: ap.Grid):
    for agent in grid.agents:
        position = grid.positions[agent]
        is_resistant_cultivar = is_intercropped_resistant_cultivar(position)
    
        # Set cell state
        agent.initial_state = AgentState.RECOVERED if is_resistant_cultivar else AgentState.SUSCEPTIBLE
        agent.state = agent.initial_state

        # Set infectivity offset based on cultivar
        agent.infectivity_offset = resistant_infectivity_offset if is_resistant_cultivar else default_infectivity_offset
        
        # Set latency duration
        agent.latency_counter = 0
        agent.latency_duration = random.gauss(default_latency_duration, 0.5)

Agent

Agent a state-machine

In the cellular automaton, each agent of the model functions as a state-machine that updates or maintains its state while the automaton's generation advances. As explained in the previous section, this is achieved by progressively calling the update_agent function to calculate the next state of every agent in the grid, based on its current state and its surrounding neighbors.

Putting the state machine in context, each vine (agent) can belong to a compartment shown in the previous section, and can transition into another compartment as it get infected by powdery mildew.

def update_agent(vine: ap.Agent, neighbors: ap.AgentList):
    if vine.state == AgentState.SUSCEPTIBLE or vine.state == AgentState.RECOVERED:
        vine.state_next = AgentState.LATENT if s_2_l(vine, neighbors) else vine.state
    
    if vine.state == AgentState.LATENT:
        vine.latency_counter += 1
        vine.state_next = AgentState.INFECTED if l_2_i(vine, neighbors) else vine.state

    if vine.state == AgentState.INFECTED:
        if self_recovery(vine, neighbors):
            vine.initial_state = vine.state_next = AgentState.RECOVERED
            vine.infectivity_offset = resistant_infectivity_offset
            return
        
        if manual_removal(vine, neighbors):
            vine.state_next = vine.initial_state
            return
        
        vine.state_next = vine.state

State transitions

This section describes the process of each state transition of the aforementioned state machine. The

Infection from neighborhood

default_infectivity = 10

def s_2_l(vine: ap.Agent, neighbors: ap.AgentList):
    neighborhood_infectivity = default_infectivity * len(neighbors.select(neighbors.state == AgentState.INFECTED))
    prob = random.random()
    return prob <= (neighborhood_infectivity / (neighborhood_infectivity + vine.infectivity_offset))

Appearing infected from latency

def l_2_i(vine: ap.Agent, neighbors: ap.AgentIter):
    return vine.latency_counter >= vine.latency_duration

Removal of infected vine

removal_threshold = 10
recovery_rate = 0.001

# TODO: make more sense out of this
def self_recovery(vine: ap.Agent, neighbors: ap.AgentIter):
    return random.random() <= recovery_rate

def manual_removal(vine: ap.Agent, neighbors: ap.AgentIter):
    return model_generation > 0 and model_generation % removal_threshold == 0

Simulation & Visualisation

parameters = {
    'size': 100,
    'initial_infection_count': 2,
    'steps': 200,
}

# Create single-run animation with custom colors
def animation_plot(model:VineInfectionModel, ax):
    attr_grid = model.vineyard.attr_grid('state')
    color_dict = {
        0:'#7FC97F', # Susceptible
        1:'#e6bc2c', # Latent
        2:'#d62c2c', # Infected
        3:'#1FC97F', # Reovered / Resistant
        None:'#d5e5d5'
    }
    ap.gridplot(attr_grid, ax=ax, color_dict=color_dict, convert=True)
    ax.set_title(f"Simulation of a SIR model\n"
                 f"Time-step: {model.t}, infected vines: "
                 f"{len(model.agents.select(model.agents.state == AgentState.INFECTED))}")

fig, ax = plt.subplots()
model = VineInfectionModel(parameters)
animation = ap.animate(model, fig, ax, animation_plot)
IPython.display.HTML(animation.to_jshtml(fps=15))

Once Loop Reflect